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A  method  of  averaging  is  used  to  study  two  general 
problems  for  the  drag  free  motion  of  an  artificial 
Earth  satellite.  First,  the  four  body  oblate  Earth 
problem,  which  describes  the  satellite  motion  per¬ 
turbed  by  the  point  mass  effects  of  the  Moon  and  Sun 
as  well  as  the  oblate  Earth  is  analyzed.  In  this 
problem  the  first  order  harmonic  J7  and  the  second- 
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order  harmonics  J,  and  J.  arc  included.  Several  trans- 
formations  are  introduced  and  the  transformed  dynami¬ 
cal  system  is  analytically  integrated.  Secondly,  two 
specific  classes  of  orbits  for  the  four  body  oblate 
Earth  problem  with  resonances  due  to  the  sectoral  and 
tesseral  potential  of  the  Earth  arc  examined.  The 
transformed  dynamical  system  for  this  preblcm  is  also 
analytically  integrated  except  for  a  numerical  evalua¬ 
tion  of  the  main  r*  sonancc  effect.  Both  problems  con¬ 
sider  satellites  with  orhital  periods  up  to  24  hours. 
The  resultant  model  is  a  mostly  analytic,  computation¬ 
ally  efficient  ephcmcris  generator  which  includes  the 
most  significant  perturbations  for  satellites  with 
orbital  periods  up  to  one  day. 


1.  INTRODUCTION. 


Four-body  oblate  Earth  treatments  include  singly 
averaged  numerical  solutions  by  Kozai  [14],  Kaula  [11]  , 
and  Cefola  [3]  ,  singly  averaged  analytic  solutions  by 
Estes  [4]  and  Koskella  and  Arsenault  [13],  and  doubly 
averaged  numerical  solutions  by  Ash  [1]  and  Lidov  [16]. 
Doubly  averaged  analytic  solutions  have  not  been  widely 
used  despite  the  apparent  potential  for  computational 
efficiency . 

There  are  very  few  references  which  treat  the  com¬ 
bined  problem  of  resonating  satellites  in  the  four  body 
oblate  Earth  problem.  The  work  of  Bowman  [2]  develops 
a  singly  averaged  numerical  model  for  12-h  elliptic 
resonating  satellites  in  the  four  body  oblate  Earth 
problem  with  drag.  A  series  of  publications  by  Kamel 
and  Tibbits  [8],  Kamel,  Eckman,  and  Tibbits  [9],  and 
Kamel  [10]  studies  a  singly  averaged,  numerically  in¬ 
tegrated  model  for  24-h  resonating  satellites  in  the 
four  body  oblate  Earth  system  by  augmenting  the  solu¬ 
tion  of  Kozai  [14].  It  is  important  to  note  here  that 
a  completely  analytic  solution  was  obtained  by  Roman- 
owicz  [20,21]  for  arbitrary  resonating  satellites  in 
the  oblate  Earth  problem.  The  Rowanowicz  development, 
however,  does  not  address  third  body  effects. 


Large  scale  space  applications  and  simulations  re¬ 
quire  efficient  and  compact  ephemeris  generators  which 
model  the  major  perturbative  effects  for  artificial 
Earth  satellites.  This  study  is  a  development  of  such 
a  generator  for  two  related  problems  for  these  satellites. 
The  first  development  applies  multiple  transformations 
to  the  four  body  oblate  Earth  problem  to  obtain  an  an¬ 
alytic  solution.  The  first-order  short  periodic  varia¬ 
tions  due  to  and  the  second-order  long  periodic 
variations  (28  day  and  365  day)  due  to  lunar  and  solar 
effects  are  recovered  using  the  transformations  intro¬ 
duced  by  a  method  of  averaging  [19]  .  The  second 
development  addresses  the  problem  of  resonating  satel¬ 
lites  in  the  four  body  oblate  Earth  problem.  These 
satellites  are  assumed  to  have  12-h  or  24-h  orbital 
periods.  The  method  of  averaging  is  again  used  to  de¬ 
rive  a  multiply  transformed  dynamical  system.  The 
transformed  dynamical  system  is  then  integrated  analy¬ 
tically  except  for  a  numerical  evaluation  of  the  main 
resonance  effect.  The  resultant  ephemeris  generators 
are  computationally  efficient.  The  second  model  re¬ 
duces  to  the  first  one  under  non-resonant  conditions. 

For  brevity,  this  discussion  is  presented  in  terms  of 
the  gravitational  functions,  and  the  explicit  expres¬ 
sions  of  the  final  solution  are  to  be  found  in  the 
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Appendices.  This  is  a  first  analysis  of  the  resonant 
four  body  oblate  Earth  problem  where  multiple  trans¬ 
formations  are  developed. 

The  result  is  a  model  which  addresses  the  same  prob¬ 
lems  as  the  referenced  papers,  but  does  not  readily 
reduce  to  any  particular  work. 

2.  THE  FOUR  BODY  POTENTIAL 

The  four  body  oblate  Earth  problem  is  a  function  of 
the  perturbative  effects  due  to  the  potentials  of  the 
Earth,  the  Moon,  and  the  Sun. 

The  spherical  harmonic  potential  of  the  Earth  is  given 
by  Kaula  [12]  as 


*  =  1 

*  sin  »*) 

where  P^is  the  associated  Legendre  functiont.  The  tonal 
potential  (R^)  can  be  separated  from  the  potential  due  to 
sectoral  and  tesseral  harmonics  (Ry)  by 

V  *  R  ♦  R 
z  v 


T+T  IW5in  5)  (Cl>n  cos 


n*0 


+  The  mathematical  symbols  are  defined  in  the  table  of 
notat ions . 
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R2  =7  a  -  J2  CtJ  P20  (5in  -  J3  (t3  p30  (sin  ^ 


(1) 


and  Ry 


a  ^ 

(-fO  P40  (sin  C)  ♦  . . .) 
is  expressed  in  terms  of  classical 


elements  by 


where  £,m,p,q  are  integers 


00 

G£pq  ^  S£mpq 

q*  -00 


i,mpq 


£-m  even 

cos  [  (ic-  2 p ) u)  +  (£-2p-*-q)M 
£-m  odd 


♦  m  (fl-9)] 


£-m  even 

sin  [(£-2p)w  +  (£-2p+q)M 
£-m  odd 

♦  m  ( fi  -  G  ^  ] 

and  (i)  and  (e)  are  functions  of  inclination  and 

eccentricity,  respectively  (see  Kaula  [12]). 

The  gravitational  potential  of  the  Moon  or  Sun  is  repre¬ 
sented  by  the  expanded  third-body  potential: 


(2) 


Rx  '  Gmx  7  '  S  p20  (cos  ♦>  ’  T  p30  (cos  ♦  >  * 

L  x  rx  rx  J 

where 

cos  $  *  cos  u  cos  cos 

-  sin  u  cos  U  cos  i  sin  (f2-H  ) 

A  A 

♦  cos  u  sin  U  cos  1  sin  (fl-H  ) 

X  X  X 

+  sin  u  sin  U  cos  i  cos  1  cos  (fi-H  ) 

X  X  X 

*  sin  u  sin  U  sin  i  sin  1 

X  X 

which  can  be  notationally  simplified  as  (see  Appendix  H) 
cos  ^  *  A  cos  f  +  B  sin  f 

A  ■  Xj  cos  F  ♦  sin 
B  *  X,  cos  Fx  ♦  X4  sin  Fx 

It  is  assumed  here  that  zonal  variations  due  to  J.  are 
first-order  effects,  while  all  other  perturbations  are 
second-order  effects.  The  24-h  limit  assumed  for  satel¬ 
lite  orbital  period  permits  a  truncation  of  the  lunar 
and  solar  potentials: 


(6) 


Rx  =  G"x  [  F  *  h  P20  (cos  *>] 

L  x  rx  J 

while  the  zonal  potential  R2  is  truncated  after  J^.  The 
potential  for  the  four-body  oblate  Earth  problem  is  then 
represented  by: 


R  =  Rz  +  Rv  +  RL  + 


R. 


The  associated  perturbed  dynamical  system  is  defined  by 
partial  differentiation  of  R  and  substitution  in  the  Lagrange 
planetary  equations  (see  Appendix  B) . 


3.  THE  EQUATIONS  OF  MOTION. 

The  complete  dynamical  system  includes  the  orbital 
motions  of  the  Moon  and  Sun  as  well  as  the  rotational  rate 
of  the  Earth,  In  this  study,  only  the  mean  anomalies  of  the 
Moon  (yl)  and  Sun  (ys)  and  Greenwich  sidereal  time  variable 

(0)  arc  assumed  to  have  significant  rates  and  it  is  assumed 

■  •  •  •  • 

that  I X"H  “GX»E  mAx»0 .  Thus  the  dynamical  system  takes  the  form 


a  3  a2  *  aL  ♦  a,.  ♦  ay 


e  »  e  *  c.  ♦  e  ♦  e 

i  L  b  V 


(i)  *  (i),  ♦  (i),  ♦  (i)c  ♦  (i), 


g  *  n2  ♦  \  -  nv 


(?) 
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(7) 


M  =  n  +  M  +M.+M  +  M 

z  L  s  v 


6  =  constant 

where  the  initial  conditons  are  defined  by  a  ,  e  . e  • 

'  O  O  r  o 

4.  THE  FOUR  BODY  OBLATE  EARTH  PROBLEM, 

a.  The  First  Transformation.  For  non^resonant  problems 
the  dynamical  system  (7)  is  a  function  of  two  fast  variables 
M  and  9,  which  satisfy  the  conditions  for  application  of  the 
method  of  averaging.  Through  Leibnitz's  Rule,  the  potential 
R  may  be  averaged  directly; 

Define: 

2tt 

<  S  >  =  yr  f  gdx 

X  4TT  J 

0 

Then 

‘  E  ‘'M.O  ■  <  RZ  >  H  *  <  RL  >  M  *  <  Rs  >  M  *  ‘  *V  >  M.8 
where 

<  "v  >  M.O  ‘  0 


and  from  Liu  ( 1 7  j  ; 


3J  .n 
4 


64a  p 


7~7  ^  ♦  jeZ)  (8  -  40  sin2 i  +  35  sin4i) 


3J3n  2 

- f-  e  sin  i  (5  sin  1  -  4)  sin  w 

8a6 


15J4n  2.2.  2  C8) 

- j-- w  e  sin  1  (6  -  7  sin  i)  cos  2u 

64a  S 


♦  3J2  "  (40  cos4i  *  8  cos2 i  -  e2  (5  -  18  cos2i 

1 2 8a20  ' 


♦  5  cos4 i )  ♦  46  (1  -  3  cos2i) 2 

♦  2e2  (1-16  cos2i  ♦  15  cos4i)  cos  2u]  , 

and 

Gill  a  **  j  y 

<  Rx  >M  -  %  ( 3A~  ♦  3B  -  2  ♦  ew  (12AZ  -  3BZ  -  3)) 

4r~ 

Tirst-ordcr  short  periodic  variations  due  to  J,  arc  given 
by  Liu  (Ref  17)  and  repeated  in  Appendix  F,  while  the 
second-order  short  periodic  variations  are  neglected. 

b.  The  Second  and  Third  Transformation.  The  averaged 
dynamical  system,  which  is  a  function  of  the  singly  averaged 
potential  <R>M  0,  is  dependent  upon  the  fast  angular 


8 


1 


variable  y^.  The  method  of  averaging  is  again  applied 
by  transforming  the  potential  to  remove  the  mean  anomaly 
of  the  Moon: 


<<  R  >vj 

M,0  yl 


<  R  >w  +  <<:  +  <  R  >.. 

z  M  L  M  y^  sM 


(9) 


The  doubly  averaged  dynamical  system,  due  to  Eq.  (9) 
is  now  a  function  of  three  angle  variables  ft,  ys,  and 
fi},  all  of  which  are  nearly  the  same  order.  The  method  of 
averaging  can  be  used  to  derive  a  transformed  dynamical 
system  which  is  independent  of  ft,  y  ,  ana  u.  This 
process  would,  however,  give  rise  to  several  singularities 
in  the  associated  transformations.  Each  singularity  would 
require  a  special  solution  in  some  small  neighborhood  of 
the  singularity.  Such  a  piecewise  treatment  results  in 
multiple  sets  of  solutions  and  is  not  compact.  Conse¬ 
quently,  this  development  uses  the  method  of  averaging 
to  derive  a  transformed  dynamical  system  which  is  indepen¬ 
dent  of  ys  alone.  The  method  of  averaging  is  applied  tc 

the  potential  -R>M  to  remove  the  periodic  dependency 

L 

upon  yt. : 


<<<  R  0>  > 
M,e  yl  y, 


<  R  >..  ♦  <<  R,  >„>  ♦  <<  R  >..> 

zM  LMyl  sMy 


S 


do) 
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where 


«  R  >M>, 


x  My, 


Gv 

8  A3  (1-Ex2)3/2  1  ~3 


(7.  +  lx  -  (4  +  6e4)) 


and 


Zx  =  3  (X2  *  X2  -  e2  (4X2  -  X2)) 


Z3  =  3  (X2  *  X2  +  e2  (4X2  -  X2)) 


2  4- 

The  long  periodic  variations  in  y^  and  yg  are  define 
by  the  method  of  averaging**  and  denoted  by: 


6R  =  — w 
x  ,  ,3 


Gmxa 


"W1  •  Ex> 


2,3/2  ^Z1F1  *  Z2F2  *  Z3F3  +  *V 


(ID 


where 


y  =  y  *  N  (t-t  ) 
1  x  x  x  o 
o 


F  = 


Yx  '  2  F-x  sin  Yx 


1  7 


*  sin  Fx  cos  Fx  *  y  (F 


X  X 


y  )  +  i  Ev  sin  F  cos2  F 
'  3  x  x  ) 


*♦  See  Ancillary  Topics  for  further  discussion, 
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ytuiaiiiiuijlltfwaiiUi 


\  E  sin  F 
3  x  x 


M  =  n  +  M  +M.+M 
z  L  s 


0  -  constant 

whore  the  secular  rates  are  presented  explicitly  in 
Appendix  D  for  zonal  expressions  and  Appendix  E  for 
third  body  effects. 

This  system  is  a  function  of  variables  which  vary 
slowly  with  time.  It  is  assumed  that  the  variable 
rates  in  system  (12)  can  be  treated  as  epoch  quantities 
for  integration  purposes,  and  that  the  resultant  integra¬ 
tion  may  be  valid  for  an  extended  period  of  time.  For 
example : 


•  •  • 


%  *  (ez  *  eL  *  es)0  (t'V 


i  “  i0  *  «*>z  *  +  ^Vo  (t’to) 


•  •  » 


■  no  *  (nz  *  \  *  Vo  '‘-'o’ 
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»  I  I 


“  *  “o  +  (“z  *  "L  *  “s>o  OV 


M  -  Mo  .  (n  *  M.  *  ml  *  Ms)0  (t-t0) 

5.  RESONANCE  PROBLEMS. 

A  resonance  condition  exists  between  M  and  0  when 
there  exist  integers  i,  m,  p,  and  q  such  that 

(£-2p)w  +  (£ - 2p+q)M  +  m(ft  -  0)  <  e  1 

where  c  < 

Equation  (13)  naturally  motivates  the  definition  of  two 


new  variables  Xj  and  X 2,  where: 


X^=M+2(fi-0)  when  M  =  20 


X^  2  u  +  M  A  ■  0  when  M  -  0 

a.  The  First  Transformation.  The  dynamical  system 
(7)  is  a  function  of  two  fast  variables  M,  and  0,  but  does 
not  satisfy  the  conditions  for  a  method  of  averaging 

■  •  i  • 

under  the  resonance  conditions  M  =  2  0  or  M  -  0.  A 
change  of  variable  is  suggested: 


M  *  X^  -  20  +  20  or  M  *  X^  -  H  -  u>  +  0 


~zrrj. sn  ~^sstz- 


which  results  in  a  new  dynamical  system, 


•  •  •  *  • 

a  =  az  ♦  aL  ♦  ag  ♦  ay 


•  •  •  •  • 

e  =e  +  e,  +  e  +  e 

z  L  s  v 


•  • 


(i)  =  (i)z  *  (i)L  *  (i)$  +  (i)y 


•  • 


Q  =  Q  +  fl.  +  +  fl 

z  L  s  v 


•  ■  • 


(l)  =CJ  +(l),  +  CJ  +  (i) 

z  L  s  v 


•  • 


X  =  Xj  orX  =^2 


VL  '  Nl 


y  »  N 
'  s  s 


(15) 


G  *  constant 

The  dynamical  system  (15)  is  now  a  function  of  one  fast 
variable,  0,  and  is  in  a  form  which  satisfies  the  condi- 
tions  of  the  method  cf  averaging.  Again  the  potential  is 
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averaged  directly; 
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where 


x  Jarctan  'WV  (*'">  evcn 

lnl  t arctan  (ctm/'stm)  (t -m)  odd 

* 

In  either  case,  Ry(X)  is  summed  over  all  possible  sets  of 
subscripts  (?,,m,p,q)  such  that  Eq.  (13)  is  satisfied. 
Therefore 

P**(X,e)  *  Ry(X,e)  -  Ry  (X)  and 

<  Ry  (X,e)  >0  =  0  and  <  Ry  (X)  >0  =  R*  (X) 
and 


<  R(x,e)  >, 


<  R  >w 
z  M 


<  RL  >M  *  <  Rs  >M  *  Rv  <X’- 


0  ’  z  M  '  "L  M  "s  M  "v  v',y  '  (17) 

is  the  singly  averaged  four-body  oblate  Earth  potential 
under  resonance  conditions.  First-order  periodic  varia¬ 
tions  in  mean  anomaly  are  recovered  as  before  while  other 
short  periodic  effects  are  neglected. 

b.  The  Second  and  Third  Transformations.  Again,  the 
singly  averaged  dynamical  system  is  in  a  form  which  per¬ 
mits  averaging  over  and  yg  successively.  The  resultant 
resonating  triply  averaged  potential  is  given  by: 


<<<  R(X,0)  >„  >  >  ■  <  R  >..>  ♦  <<  R,  >u> 

9  yL  ys  z  M  1  M  yL 
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_ _ 

.  -  ^  *■  fm? 


1 


+  <<  R  >M>  +  R  (x) 

S  M  y  v  v  1 

and  the  associated  dynamical  system  has  the  form: 


a  =  a 


•  • 


e  =  e  +  e,  +  e 

z  L  s 


(i)  =  U)z  ♦  (i)L  +  (i) 


0  =  +  h  +  ns 


U  =  +  WL  +  us 


•  • 


X  ■  X,  or  X  -  X2 


yl  ■  nl 


Y  *  N 
1  s  s 


(18) 


0  M  constant 

"  *  *  *  *  *  *  *  ’  * 

where  ey,  (i)v,  ny,  uv,  and  My  are  neglected.  The  variable 

*  *  * 

a  may  be  obtained  by  partial  differentiation  of  Ry  (X)  and 

substitution  into  the  Lagrange  planetary  equation  for  a. 
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It  is  convenient  at  this  time  to  change  variables  and  re- 
place  a  by 

n*  =  (-?n/2a)av 

c.  Four-Body  Integration.  The  dynamical  system  (18) 
is  again  a  function  of  variables  which  vary  slowly  with 
time,  except  n  and  A,  which  may  oscillate  rapidly.  Ex¬ 
cept  for  n  and  X,  the  orbital  element  rates  in  the  sys¬ 
tem  (18)  may  be  treated  as  epoch  quantities  for  integration 
purposes.  Further,  most  of  tiie  variable  expression  in  n 
and  A  are  assumed  to  be  epoch  constants,  except  for  n,  w, 
and  X  as  follows: 

nv  =  (  3noao/w)  S  (QimpqJo  Sin  (cu> 

i »m,p,q 


♦  (fc - 2p+q)  X  -  A,  ) 

(19) 

•  •  •  •  •  •  • 

n  *  (Hj  *  Ml  •  Ms  *  '(B,  *  BL  *  ns)-26)0 

««•••••••• 

"  *  (M.  *  ml  *  Ms  *  q2  *  at  ♦  ns  ♦  -  m5  -  s)0 


where  c  a-q  or  ( S- - 2 p)  (see  Eq.  (16))  and  u  *  +  u>q  (t-tQ) 

for  efficiency,  a  numerical  integration  of  Eq.  (19)  is 
chosen,  resulting  in  a  mostly  analytic  integration  cf 
the  system  (18).  It  should  be  noted  that  the  argument 
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of  perigee  must  be  allowed  to  vary  with  time  during  the 
integration  of  Eq.  (19),  particularly  for  the  highly 
eccentric  12-h  satellites. 


6.  ANCILLARY  TOPICS. 

a.  Periodic  Expressions.  The  periodic  variations  due 
to  and  ys  are  recovered  in  the  discussions  of  second 
and  third  transformation.  These  periodic  variations  are 
functions  of  the  integrated  transformed  variables  and  are 
quite  complex.  To  be  consistent  with  the  assumptions 
made  in  the  integration  sections,  all  variables  except 
and  Fs  are  evaluated  under  epoch  conditions.  Thus: 
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Both  long  and  short  periodic  variations  arc  singularity 
free  through  a  change  of  variables.  An  example  of  such 
a  change  of  variables  is  to  be  found  in  Appendix  G. 

b.  Small  Divisors.  The  classical  element  formula¬ 
tion  of  the  secular  systems  (12)  and  (18)  suffers  singu¬ 
larities  for  zeros  in  inclination  and  eccentricity  in 
second  order  expressions  from  R.  ,  R_,  and  R^  (J-  only), 
these  systems  can  be  formulated  in  terms  of  equinoctial 
elements  (sec  Ref  (3)).  The  expressions  with  classical 
inclination  singularities  hecomc  third  order  non-singular 
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expressions  in  a  neighborhood  of  the  singular  point. 

The  classical  eccentricity  singularities  occur  in  the 
terms.  For  high  altitude  satellites  where  lunar  and 
solar  perturbations  are  second-order,  the  perturbations 
due  to  J,  become  third-order  when  eccentricity  is  small. 

In  either  case  the  expressions  with  classical  singularities 
become  third-order  near  the  singular  points  and  can  be 
neglected  in  this  second-order  development.  The  same 
result  is  more  economically  obtained  by  simply  neglecting 
the  singular  classical  variable  expressions  in  some 
small  neighborhood  of  the  singularity. 

c.  The  Ephemeris  of  the  Moon  and  Sun.  The  classical 
orbital  elements  of  the  Moon  and  Sun  with  respect  to  the 
equatorial  plane  at  the  time  of  epoch  of  the  satellite 
elements  are  obtained  using  the  model  presented  in  the 
Explanatory  Supplement  to  the  Astronomical  Ephemeris  and 
American  Ephemeris  and  Nautical  Almanac  15].  The  model 
is  supplied,  wit'.. out  explanation  in  Appendix  I. 

d.  Other  Simplicat ions .  Through  order  of  magnitude 
analysis  of  the  Qimpq  coefficients  in  Eq.  (19),  a  signi¬ 
ficant  reduction  in  the  number  of  terms  in  the  resonating 
potential  is  possible  depending  upon  the  desired  accuracy. 

For  the  existing  2 4 - h  satellites,  the  most  significant 
^tmpq  arc  defined  by  the  set  of  quadruplets  (£,m,p,q)  ■ 

{(2,2, 0,0),  (3, 1,1,0),  and  (3, 3, 0,0)}  .  For  the  existing 
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1 2 - h  satellites,  the  set  of  significant  quadruplets  is 
much  larger  due  to  the  ranges  of  values  for  inclination 
and  eccentricity.  In  general,  the  more  significant  12-h 
quadruplets  arc  those  where  i  <  6  and  m  <  6  and  Eq.  (13) 
is  satisfied.  The  results  obtained  for  this  paper  are 
based  upon  these  simplifications.  For  NO RAD  applications 
the  generalized  G^q(e)  and  ( i )  •  as  presented  in 

Appendix  J,  proved  to  be  too  cumbersome.  Explicit  ex¬ 
pressions  have  been  developed  for  the  most  frequently 
used  eccentricity  and  inclination  functions;  these 
explicit  functions  are  also  presented  in  Appendix  . 

7.  RESULTS. 

The  graphical  results  presented  here  are  long-term 
comparisons  versus  epoch  orbital  elements  of  existing 
satell  tes  extracted  from  the  NORAD  Historical  Data 
System.  (Data  anomalies  due  to  historical  circumstances 
are  not  edited).  The  values  of  the  N0RA1)  elements  are 
represented  by  X's.  The  solid  lines  arc  predictions 
with  this  ephemeris  generator,  assuming  the  "X"  denoted 
by  "H"  as  the  cpo:h  value.  The  SAO  1969  geopotential 
model  is  used. 

Figures  1  and  2  depict  the  mean  motion  and  east 
longitude  (1,)  for  the  24-h  satellite  1971-095B  for  the 
period  September  1972  through  April  197S.  Resonating 
about  both  "null  points",  this  atellitc  almost  cir¬ 
cumnavigates  the  Earth. 
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Figures  3  and  4  present  the  mean  motion  and  resonance 
angle  (Xj)  for  the  elliptic  12-h  satellite  1964-049E 
for  the  period  January  1970  through  July  1971.  This 
satellite  represents  an  extreme  case  of  resonating  12-h 
satellites  in  terms  of  maximum  drift  rate  in  Xj  (approxi¬ 
mately  l?35/day) . 

These  two  test  cases  provide  a  successful  examination 
of  this  model  for  resonating  satellites  in  the  four-body 
oblate  Earth  system. 

8.  REMARKS. 

Throughout  the  region  of  interest  it  is  assumed  that 
the  zonal  effects  due  to  J^  are  first-order  variations 
while  all  other  perturbations  are  second-order  variations 
This,  of  course,  may  not  be  valid  for  cases  in  which 
the  lunar  and  solar  effects  are  either  larger  or  smaller 
than  second-order. 

For  near-Earth  satellites,  where  lunar  and  solar  per¬ 
turbations  are  third-order,  the  reader  is  referred  to 
one  of  the  many  oblate  Earth  second-order  theories. 
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10.  NOTATION. 


All  variables  are  referenced  to  a  non-rotating 
Earth-centered  coordinate  system  defined  by  the  equatorial 
plane  and  the  line  of  poles.  In  this  system,  the  Sun 
"appears"  to  revolve  about  the  Earth  and  is  thus  given 
orbital  elements. 
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I  . 

1 

F 

x 

F 

imp 

f 

G 


Definition 

Intermediate  variable 
Semimajor  axis  of  body  X 

Semimajor  axis  of  satellite 

Mean  equatorial  radius  of  the  Earth 

Intermediate  variable 


Cosine  coefficient  of  spherical  harmonic 
potential 

Eccentricity  for  body  X 

Eccentricity  of  satellite 

Periodic  functions  of  Fx  (see  Eq.  11) 

True  anomaly  of  body  X 

Inclination  function 

True  anomaly  of  satellite 
Gravitational  constant 
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Symbol 

g 

H 

x 

I 


J  . 
1 


L 

l 

M 

m 

m 

) 

N 


fcm 


£mpq 


q 

R 


r 

r 

x 


Definition 
Arbitrary  function 
Longitude  of  node  of  body  X 

Inclination  of  body  X 

Inclination  of  satellite 
Coefficient  of  zonal  potential  term 

Subscript,  referencing  the  lunar  potential 
Subscript,  degree  of  spherical  harmonic 
Mean  anomaly  of  satellite 
Subscript,  order  of  spherical  harmonic 
Mass  of  body  X 

Mean  motion  of  body  X 

ilean  motion  of  satellite 
Subscript,  evaluated  at  epoch 
Legendre  associated  polynomial 

Subscript,  index  of  inclination  function 
Amplitude  of  term  in  spherical  harmonic 
potential  disturbing  function 
Subscript,  index  for  eccentricity  function 
Disturbing  function 

Periodic  expression  from  method  of  averaging 

Radial  distance  to  satellite 
Radial  distance  to  body  X 
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Symbol 

s 


Jlmpq 


t 

U 

u 

V 


X 

Z. 

1 

6 

e 


v 

C 

Y 

$ 

n 

u 


Definition 

Subscript,  referencing  the  solar 
potential 

Sine  coefficient  of  spherical  harmonic 
potential  term 

Combined  spherical  harmonic  potential 
term 

Time 

F  +  Gv 
x  x 

f  +  (i) 

Subscript,  referencing  the  Earth's 
harmonic  potential 

Intermediate  variable 

Subscript,  read  L  for  lunar;  s  for  solar 
Intermediate  variable 

^  1  -  e2 

Greenwich  sidereal  time 
Resonance  angle(s) 

Gravitational  constant  of  Earth's  mass 
Latitude 

Mean  anomc'y  of  body  X 
Phase  angle  between  r  and  r 

Longitude 

Longitude  of  node  of  satellite 
Argument  of  perigee  of  satellite 
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APPENDIX  A 


A  REVIEW  OF  THE  METHOD  OF  AVERAGING 


The  following  discussion  is  extracted  verbatim  from  Liu 
[18]  with  permission  of  the  author  and  formally  documents 
the  procedures  of  the  Method  of  Averaging.  For  a  more  com¬ 
plete  discussion,  readers  are  referred  to  Morrison  [19]  or 
Kyner  [ 15] . 


Consider  a  system  of  first-order  ordinary  differential 
equations  written  in  the  component  formt 


xi  cXu  <v 


yn)  *  e  x2i  (X 
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-  yv  *  rYu  <v  v  *  e2y2j  <v  v 

where  i  =  l,2,.,.,M,j  =  1 ,2 , . , . ,N ,  1  £  m  £  M,  and  l£n<N 
with  the  initial  conditions  xm ( 0 j  =  am,  y n ( 0 )  «  bfi.  The 
x  ,  are  referred  to  as  slow  variables  because  their  time 
variations  are  proportional  to  the  small  parameter  c.  The 
y  arc  referred  to  as  fast  variables  because  the  dominant 
parts  of  their  time  variations  are  proportional  to  t.  The 
functions  Xu  Um;  y„ ) .  X2i  (x,;  y„)  ,  Y,j  lx„.  >„)  . 


+  If  F . ,  i  -  1,2,...L,  are  functions  of  the  variables  x^, 

m  *  1,2,...,M,  we  write  F^(xm).  If  F^  arc  functions  of 

variables  x  and  y  .  m  *  l,2,...,M,n  *  1,2,. ...N,  we  write 
m  n 

Hi  <V 

A-l 


and  Y, .  (x  ;  )  are  assumed  to  be  continuous  functions  of 

l  j  m  n 

xm  and  yn  with  a  period  of  2 tt .  In  general,  the  dynamical 
system,  Eqs.  (A-l),  is  non-linear  in  nature  and  complex  in 
form  and  hence  the  integration  of  the  system  is  usually 
analytically  intractable. 

As  the  first  step  in  the  solution,  a  transformation  is 
introduced 


; .  =  x.  +  eP..  (x  ;  y)  +  c  ZPV  (x  ;  y  ) 
i  1  li  v  m’  Jn  2i  v  m  7nJ 


y.  =  y  ♦  eQ,  •  (x  ;  y  )  +  c  Qv  (x  ;  y  ) 
})  xlj  v  m’  7nJ  N2j  v  m’  'n 
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so  that  in  a  sense  the  differential  equations  become 
simpler  to  handle  in  terms  of  the  new  variables  x  and  v  . 

Here  x.  and  y.  arc  regarded  as  the  new  unknowns  and  .  (xm;  yp) , 

P2i  <V  V'  Qlj  <V  «2j  txm;  are  pcriodic  func- 

tions  of  each  y  as  new  functions  to  be  determined  in  such 
a  way  as  to  effect  a  simplification  in  the  transformed  dynami¬ 
cal  system.  It  is  desired  that  the  fast  variables  y  (to 
the  second-order  in  c)  be  eliminated  from  the  transformed 
differential  equations.  Thus,  the  transformed  differential 
equations  are  to  have  the  form 
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For  a  second- 


for  suitable  functions  U^,  ,  ^2j* 

order  theory,  the  third-order  terms  and  will  be 
ignored.  Thus  Eqs.  (A<-3)  become 


xi  =  EUli'xJ  *  6 


=  Zj(xm>  *  *  E  V2jtxm) 


(A-4) 


with  initial  conditions  x  (0)  =  a  .  y  (0)  =  b  .  These 
initial  values  are  obtained  by  substituting  the  initial 
values  am,  bn  into  the  transformation  (A-2).  The  explicit 
expressions  of  functions  IJ^,  V^,  v2 j  ’  Pli’  P2i* 

,  and  Q2j  will  be  given  here  without  proof.  For  the 
thorough  discussions  of  the  general  theory  of  the  method  of 
averaging,  see  Refs  19  and  15. 

To  describe  the  explicit  expressions  of  these  functions, 
the  necessary  relations  and  their  definitions  will  be 
introduced.  To  begin  with,  the  Fourier  series  expression 
for  the  perturbing  function  has  the  form 
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In  Eqs  (A-6),  the  following  definitions  apply.  The 

notation  k  =  (kj  ,  V.2 . kn . . k^)  is  used  to  denote  a 

vector,  each  component  k^  is  an  integer;  k  =  0  is  not  per* 
mitted.  The  notation  j[  indicates  summation  over  all  possible 
integer  vectors  k. 

The  notations  X2ilkc^xm^  and  Xlilks^xm^  denote  the  coefficients 

of  the  cos  [k,y]  and  sin  [k,y]  terms,  respectively,  in  the 
summation.  These  Four ier- series  coefficients  are  to  be 
determined  in  the  usual  way.  It  can  be  shown  by  the  method 
of  averaging  that 

Uli(xm'1  "  XliO(xrJ’  V1  j  *xm*  “  ’'ljO^m^ 
and 
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In  Eqs .  ( A-7)  and  (A-8),  the  notations  Y1;j0(xm)»  Y^1  (xm;  yn)  , 
Sljks^*m^  an<*  ^Ijkc^m^  ^ave  t*ie  similar  definitions  which 
apply  for  *aoljJ.  XH1(in;  y„) ,  Xj^tfJ  and  Xj^tf,). 
respectively.  Secondly,  two  new  functions  are  introduced 
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The  method  of  averaging  shows  that 


U2i(xJ  ■  RliO<xm>'  V2j(xm>  "  R2jOtxm> 
and 


(A- II) 
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Again  the  notations  Rli0(xn),  \n  (*V  y„) ,  *ni±si*m) , 

Rlilkc(xm):  R2j0^xm^»  R2jl(xm;  yn^  *  and  52jks’  Ssjkc  V 
have  the  usual  meaning  in  the  Fourier  series  expressions 
for  and  ,  respectively.  It  should  also  be  noted, 

in  arriving  at  these  expressions,  that  the  nonresonant  condi¬ 
tion,  i.e.,  (k,Z(xm)]  +  0,  has  been  assumed. 
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Following  the  determination  of  the  functions  U's,  V's, 

P's,  and  Q's,  Eqs.  (A^4)  may  either  be  integrated  analytically 
or  integrated  numerically  with  a  much  longer  integration 
step  time  than  may  be  used  with  Eqs.  (A-l).  The  second- 
order  solution  for  and  y.  can  then  be  obtained  from  Eqs. 

(A- 2) . 
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APPENDIX  B 


LAGRANGE  PLANETARY  EQUATIONS 

The  Lagrange  Planetary  Equations  are  presented  here  for 
completeness : 


a  = 


2na  3R 

T”  m 


e  = 


1_ 

ae 


P* 


1/2 


3R  np*  3R 
TO  ye  TO 


(i)  -  - 


sin 


3R  .  cot  i  9R 

T  TO  H  TO 


n  = 


_ 1 _ ^  3R 

h  sin  1  TO 


u  = 


Sie  IS  '  (cos  “ 


M 


np*  3R  _  2na^  /  3R  \ 
ye  TO  y  \  TO  / 


M 


where 


h  *  na^fi 
p*  *  aB1- 


B- 1 


APPENDIX  C 

THE  TRUNCATED  THIRD  BODY  POTENTIAL 

1.  The  truncated  potential  for  a  third  body  (x)  is  given 
by 

2 

Rx  =  Gmx  P2  ^cos  ^ 
rx 

where 

cos  4>  =  cos  u  cos  U  cos  (ft  -  H  ) 

X  X 

-  sin  u  cos  i  cos  U  sin  (ft-  H  ) 

X  X 

+  cos  u  sin  U  cos  I  sin  (ft  -  H  1 
x  x  x 


+  sin  u  cos  i  sin  U  cos  I  cos  (ft  -  H  ) 

x  x  v  xJ 


+  siri  u  sin  l  sin  U  sin  I 

X  x 


C-l 


2.  The  averaged  potential  for  a  third  body  is  given  by: 

Rx  =  — — — —  [3A2  +  3B2  -  2  +  e2  (12A2  -  3B2  -  3)] 

4rx 

where : 

A  =  cos  w  cos  U  cos  (ft  -  H  ) 

X  X 

-  sin  u>  cos  i  cos  sin  (ft  -  Hx) 

+  cos  w  sin  U  cos  I  sin  (ft  -  H  ) 

X  X  v  X 

+  sin  oj  cos  i  sin  U  cos  I  cos  (ft  -  H  ) 

X  x  v  XJ 

+  sin  w  sin  i  sin  U  sin  I 

X  X 

B  *  -  sin  td  cos  Ux  cos  (ft  -  Hx) 

-  cos  <d  cos  i  cos  U  sin  (ft  -  H  ) 

DC  DC 

-  sin  (d  sin  U  cos  I  sin  (ft  -  H  ) 

X  X  v  X 

+  cos  (d  cos  i  sin  U  cos  I  cos  (ft  -  II  ) 

X  X  X 

♦  cos  u)  sin  i  sin  Ux  sin  Ix 

3.  The  potential  for  a  third  body,  averaged  first  over  one 
period  of  the  satellite,  and  further  averaged  over  one 
apparent  period  of  the  third  body  is  given  by: 


i 

i 

} 


C-2 


where  Z^,  are  given  in  Appendix  H, 

2 

(note  that  ( 1 - Ex  )  terms  are  dropped  in  NORAD  applications) 


APPENDIX  D 


ZONAL  SECULAR  EFFECTS 

The  singly  transformed  dynamical  system  derived  by 
Liu  [17,18]  is  given  in  terms  of  classical  elements. 

a_  =  0 


ez  =  -  n J 2 2  (  p  j  sin^i  (14-15sin2i)  e  (1-e2)  sin  2ii) 


-  ^  nJ^  (  ^  )  sin  i  (4-5sin2i)  (1-e2)  cos  u> 


15  nj  /  R  sin2i  (6-7sin2i)  e  ( 1  - e ^ )  sin  2  u> 
37  J4  l  p 


3  T  2  /  R  sin  2i  (14-15sin2i)  e2  sin  2  u 
[1)z  =  M  nJ2  Ip 


/  \  3 

♦  nJ3  |  5.  |  cos  i  (4-5sin‘‘i)  e  cos  w 


2 . 


+  ^  nJ4  l  ~  ^  sin  2i  (6-7sin2i)  e2  sin  2  u> 


2, ,  .  3  _  T  2  /  R 


-z  -  TnJ2?  C4’5sin  iJ  +  KnJ2  U  * 


2  21  ^  4  /  9  2  45.4. 

48-103sin  l  ♦  ^  sin  i  ♦  I  7  *  y  sin£l  -  ^  sin  1 


D-l 


'  *  «•'  >■*-  ---i  J  E  — *«■ 


+  611 


^1  *■  j  sm  ij  (4-5sin  i)  (l^e  )  ' 


1  [2(14-15sin2i)  sin2i  -  (28^158sin2i 


+  135sin^i)  e2]  x 


3  /  R 

cos  2  to  •+  it  nJ_  — 
8  3  Ip 


(4-SsinZi)  ^nZi-e2cos2i 
e  sin  i 


+  2  sin  i  (13-15sin2i)  e 


j  sin  a)  - 


if  nJ4  (  |  )  x 


|  16-62sin2i  +  49sin^i  +  (24-84sin2i  +  63sin^i)  e2 


+  sin2i  (6-7sin2i)  -  ^  (12-70sin2i 


+  63sin^i)  e2  cos  2  a) 


nJz{l  )  cos  1  -  1  nJ2Z(f  )  cos  1  If 


*  |  (l-c2)1/Z 


-  -g-  1  +  5sin2i  -  sin4ij  e4  +  sin2i  x 


[70- 123sin2i  +  (S6-66sin2i)  e2]  cos  2  u 


+  TIS  s^n4i  cos  4  w  j  -  g-  n^3  (  p  ^  s^n  4  (4-5sin2i) 


l-4e  ,,  2,1/2  45  T  /  R  \ 

—r~  (1‘e  ]  sin  w  '  m  nJ.<  u  (8 


4Csin2 i 


+  35sin4i)  e2 


^  -SMS) 


sin2i  (6-7sin2i)  (2-5e2)  (1-e2)^2  cos  2  u 


where  R  is  the  equatorial  radius  of  the  Earth,  and  p  -  a(l-e  ) 


APPENDIX  E 


THIRD  BODY  SECULAR  EFFECTS 

The  following  equations  represent  the  secular  effects 
due  to  a  third  body  (subscript  x)  .  These  expressions  can 
be  obtained  by  partial  differentiation  of  the  doubly  averaged 
third  body  potential  as  given  in  Appendix  C.  As  noted  in 
the  text,  the  right  hand  sides  are  assumed  to  be  constant 
functions  of  the  epoch  elements  of  the  satellite  and  the 
third  body.  The  equational  forms  of  the  and  Z ^  are 
given  in  Appendix  H. 


if  i  <3°  set  (D  sin  i)^  -  0 


E- 1 


ft  =  (ft  sin  i)  /sin  i 

A  A 


u  =  (to  +  ft  cos  i) 

x  '  X 


cos  i  ft 


E-2 


- . -  -  -  - 


_ _ 


APPENDIX  F 


6fl. 


e  cos  (2u+f)  +  cos  2  (tu+f)  +  |  cos  (2w*3f)j 

jf  -  M  +  e  sin  f  -  |  sin  (2u+f) 


?  "2  l  p 


COS  X 


1 


{“z  -  f  J2 


f  sin  2  (ui+f)  -  g-  sin  (2w+3f) 


(  f  )  <4  '  5sin2i)  (f  -  M  ♦  e  sin  fj 


‘  1  J2 


«)'  ( 


,  3.2. 

1  -  £  sin  i 


OfJM-’) 


sin  f 


♦  |  sin  2f  ♦  i  e  sin  3f 


sin2i 


^  1  -  g  sin2i j  sin  (2w+f)  +  sin2i  sin  (2w-f) 


5  sin2i  )  sin  2  (ut+f)  -  i 


*  \  (  1  *  f  sm2i  ^ 

’  T  (  1  '  X  sin2i  )] 


7  .2. 

u  bin  1 


sin  (2w+3f) 


3  2  12 

8  sin  i  sin  ( 2u)+4 f )  -  e  sin  i  sin  (2w*5f) 


i_  , 
16  2 


(?) 


•  2 

sin  i  sii.  2<d 


F-  2 


MMMTt 


X 


5Mz  ■  ■  I  J2  j? 


-e2l1/2  j  /,  3.2. 

— e  j  (1  '  2  Sln  1 


^  1  -  j  e2  ^  sin  f  +  |  sin  2f  +  yy  s:*-n  3f 

+  j  sin2i  |^'7(1+|'e2)  sin  (2a)+f) 

2  7  /  2  . 

-  sin  (2w-f)  +  g  (  1  -  jg-  )  sin  (2u+3f) 


+  -r  e  sin 
4 


(2u+4f )  +  ^  sin  (2<j+Sf)  j  [  +  ^  J2  (  p  )  (1 


-  e2)1^2  sin2i  sin  2w 


where  R  is  the  equatorial  radius  of  the  Earth  and  p  =  a(l-e  ) 


APPENDIX  G 

THIRD  BODY  PERIODIC  EFFECTS 


The  following  equations  represent  the  additive  Lunar 
(second)  inverse  transformation,  as  given  in  Eq,  (20)  of 
the  text,  where  X.,  Z . ,  F.  and  C.  are  as  given  in  Appendix 
H  for  Lunar-related  variables.  To  obtain  the  Solar  (third) 
inverse  transformation,  determine  the  X.,  Z.,  and  C  for 

1  1  S 

the  epoch  solar  position  and  F^  at  the  time  of  interest  and 
substitute  below: 


6eL  = 


- 30&CLeo 

n 


(F1X1X3  +  F2  (X2X3  +  XiX4 )  +  F3X2x4) 


5lL  -  n  6  ^F1Z11  *  F2Z 1 2  *  FJZu' 


{ML  =  4^  (F1Z1  *  F2Z2  *  FiZJ  *  3EL  5in  FL  ('7  '  3  eo» 
o 

2  c 

SghL  •  (FjZjj  *  F2Zj2  *  FjZjj  -  9  El  sin  FL  ) 

o 


6hL  =  n  B  (F1Z21  +  F  2 Z  2  2  +  F3Z23^ 


The  long  period  third  body  terms  for  0  and  w ,  when 
i  >  .2  radians,  are  computed  by: 

6uj^  *  <5gh^  '  cos  ^^L^  s*n  io 


G-l 


fin  a  fihT /sin  i 

L  O 


However,  as  reported  in  the  text,  the  divisor  by  sin  i 
for  small  inclinations  necessitates  the  following  change 
to  intermediate  variables: 

fia^  =  fih^  cos  ft  +  <5i^  cos  i  sin  ft 

fiB*^  =  <5hL  sin  ft  +  fii^  cos  i  cos  ft 

fiL*^  =  +  =  Si^  ^  sin  i 

where 


L*  =  Mp  +  to  +  ft'cos  (i) 


a  *  sin  (i)  sin  ft 

3*  »  sin  (i)  cos  ft 

and  where  the  variables  on  the  right  hand  sides  are  inte¬ 
grated,  third  transformed  variables. 

In  this  special  case  the  periodics  are  applied  through 
the  intermediate  variables  L* ,  a,  and  6*.  The  integrated 


G-2 


singly  transformed  variables  are  denoted  by  the  p 
subscript : 


a 


a  +  6aT  +  6a 
L  s 


0*  =  0*  +  60*  +  60* 

L  S 


L*  =  L*  +  61*  +  6L* 

L  s 

ftp  =  tan  1  (u/0*) 


i  =  i  +  6i.  +  6i 
p  L  s 


Mp  =  Mp  +  6Ml  +  6Ms 


u  =  L*  -  M_  -  ft„  cos  (i  ) 
u  p  p  .> 


e 

P 


e  +  6eT  +  6e 
L  s 


L 

P 


M  +  id  +  ft 
P  P  P 


G  -  3 


- - — — 


'-**1**.' — ■*-  *-  "  'Ur 


APPENDIX  H 

INITIAL  VARIABLES  COMMON  TO  APPENDICES  E  AND  G 

The  following  variable  calculations  are  intermediate 
and  in  themselves  have  no  important  physical  definition. 
This  collection  of  terms  was  chosen  to  minimize  the  number 
of  calculations  in  obtaining  third  body  effects. 


al  =  cos  Go  cos  (no  "  H0  5  +  sin  Go  cos  lo  5in 


-  H  ) 
o 


a,  =  -  sin  G  cos  (n  -  H  )  +  cos  G  cos  I  sin  (n 
3  o  oo  o  o  vwo 

X  XXX 


-  H  ) 
0 

X 


cos  G  sin  (o  -  H  )  +  sin  G  cos  I  cos  (n 

n  '  ■*  q  O  n  n  v 


°x  °x  0 


Ho  > 

X 


a0  =  sin  G  sin  I 
8  °x  °x 


a^  =  sin  Gq  sin  (Dq  -  Hq  )  +  cos  Gq  cos  Iq  cos  (dq 


-  Ho  J 

x 


a  •  cos  G  sin  I 

X  X 


H- 1 


a_  cos  1  +  a„  sin  i 

7  o  8  o 


a4  =  a9  cos  xo  +  a10  Sin 


ar  =  -a,  sin  l  +  an  cos  i 
5  7  o  8  o 


a,  =  -an  sin  i  +  a. _  cos  i 
69  o  10  o 


X.  =  a.  cos  a)  +  a_  sin  u 
11  o  2  o 


=  a,  cos  u)  +  a.  sin  u 
2  3  o  4  o 


XT  =  -a,  sin  w  +  a_  cos  w 
3  1  o  2  o 


X,  =  -a_  sin  ui  +  a.  cos  w 
4  3  o  4  o 


Xr  =  ar  sin  u 
5  5  o 


X,  =  a,  sin  w 
6  6  o 


X_  =  ac  cos  u 
7  5  o 


X0  =  a,  cos  u) 
oo  o 


"31  ■  12X1  '  3X32 


Z-,  *  24X, X_  -  6X_X . 
<j  ^  1  Z  j  4 


Z33  '  ,2X2  •  3X4 
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21  =  3X1  +  3X3  +  eo  231 


22  -  eX.^  ♦  6X3X4  +  232 


Z3  -  3X2  *  3X42  ♦  e2  23J 


2..  ■  -6X.X-  +  6X  X.  ♦  e  (-24X.a„ 
11  17  35  o  17 


2 12  *  ‘6XiX8  -  6X2X7  +  6X3X6  +  6X4X 


►  e„  ( - 24X_X_  -  24  X,  X0  -  6X,X, 
o'-  2  7  1  8  3  6 


2 13  =  -6X2X8  *  6X4X6  +  eo  ( ‘ 24X2X8 


22i  =  6XXXS  ♦  6X3X?  ♦  ( 24X1X, 


2,,  =  6X,XC  ♦  6X, X  +  6X.X-  +  6X,X0 
22  25  16  47  38 


+  (24X,XC  ♦  24X.X,  -  6X-X- 

o'-  2  5  1  6  4  7 


2 ,,  -  6X_X ,  +  6X.X0  *  e  (24X,X. 
2j>  26  48  o  26 


Redefine  2 j ,  Z^,  and  2^  in  terms  of  2. 


2 1  =  22i  ♦  (1  -  e‘)  23i  i  =  1,2, 


The  constants 


C,  ■  4.7968065E-7  rad/min 


Cs  =  2 . 98647969E-6  rad/min 
are  derived  from 

]_ _ x _ 

C  =  4  m  +  m  N 

A  A  g  A 


where 


m^  =  mass  of  perturbing  body  (Moon  or  Sun)  and 

m  =  mass  of  the  Earth 
e 


Nx  =  apparent  mean  motion  of  Moon  (or  Sun)  about  the  Earth 
For  the  Sun  we  choose  to  approximate 


m  +m 
x  e 


while  for  the  Moon  we  use  m  *  6 . 1  r ,  in  . 

x  81.53  e 

The  mean  motions  of  the  Sun  and  Moon  are  given  by  (see 
Appendix  1  also) : 


Ns  =  1.19459E-S  rad/min 
Nl  =  1. 583521770F.-4  rad/min 

The  variables  F^,  F^,  and  F^  are  conputed  as  follows  where 
subscript  x  *  L  or  s  for  Moon  or  Sun: 


H-4 


(The  F  presented  here  are  a  simplification  of  the 
presented  in  the  text  (Eq.  11)). 


H*S 


APPENDIX  I 

LUNAR  AND  SOLAR  EPHEMERIDES 

The  positions  of  the  Moon  and  Sun  in  terms  of 
classical  equatorial  elements  are  given  as  follows  (see 
ref  [5]): 


ftL  -  259° 1833275  -  0? 0 5 2 9 5 39 2 2 2 t * 
Le 


cos  i.  =  cos  e  cos  iT  -  sin  e  sin  iT  cos  G. 
Lo  Le  Le  Le 


•  /  ,  2. 

sin  iT  =  /  1  -  cos  i_ 

Li  Li 

o  o 


270 ?4342  ♦  13°17oi965268t* 


334?32955  ♦  0 ?1114040803t* 


Y  -  m*  -  t; 
o  o 


sin  Ht  sin  iT  -  sin  i_  sin  GT 

Li  JL  Li  Li 

o  o  c  e 


cos  H,  *  1  -  sin  H_ 


♦  tan 


sin  A 
cos  A 


Y  *  358?4 7  584  +  0?985600267t* 
So 

Nl  =  13?06499244/day 


Ng  =  ?98S600267/day 

where  t*  is  the  time  in  days  since  J.D.  241  5020.0  (1900  Jan. 
0.5),  and 

iL  -  5?1 4 5 396 
c 

e  =  23?4  4  4 1 
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APPENDIX  J 

INCLINATION  AND  ECCENTRICITY  FUNCTIONS 


The  generalized  eccentricity  and  inclination  formulae, 
G^pqte)  and  Ftmp^’  are  8iven  bclow  (see  Kaula  [12]), 

These  equations  are  extracted  verbatim  from  Kaula  for 
completeness.  For  terms  in  which  fc  -  2p  +  q  =  0 ,  compute: 


G£p(2p-1) "  (1_e2}£- 


p  y 1  (  i-1  \ 

1  I  2d  +  £  -  2p' 


2d  +  %  - 2p ' 
d 


(?) 


2d+t - 2p ' 


in  which 


p'  -  p  for  p  <  1/2, 


p'  =  £-p  for  p  >  l/l. 


For  the  terms  where  t  -  2p  q  /  0,  the  development  of 
Glpq(e)  *s  muc^  more  complicated;  wc  merely  quote  the  result 
of  one  solution  (Tisserand,  1889,  p.  256): 


Gtpq(e)  -  (-1)  (I  *  bYs""  ^  ?t[)qkQtpqk82k 


h  =  k  +  q  ’  ,  q'>0;h“k,q'<0; 


and 


h  ■  k,  q  ’  >  0;  h  *  k  -  q',  q  ‘  <  0; 

p  *  p,  q’  *  q  for  p  <  1/ 2;  p'  -  £-p,  q*  =  -q  for  p  >  1/2 
While  the  inclination  function  p£mpO)  is  given  by 


F£n,p^ 


l  _ C Z*-ZQ!  , 

t  t! ( £-t) t (£-n-2t) 12‘ 


.  I  m-2t . 
sm  l 


-  2t 


Here  k  is  the  integer  part  of  (j.  -  m)/2,  t  is  summed  from 
0  to  the  lesser  of  p  or  k,  and  c  is  summed  over  all  values 
making  the  binomial  coefficients  nonzero. 


As  discussed  in  the  text,  the  G^^Ce)  expressions  are 

much  too  cumbersome  and  too  computationally  slow  for  a  fast 

general  perturbations  theory  because  of  the  eccentricities 

with  which  we  are  dealing,  while  the  F0  (i)  have  an  exact 

mp 

formulation.  The  F^^Ci)  are  used  as  follows; 


F  2  2 1  *  I  ^Sin  ^ 

F321  *  T  sin  U0>  ^  '  2  C0S  (io)  '  3  cos2  (i0» 

F322  *  ‘  X  sin  ^o5  ^  +  2  cos  (V  ’  cos2  (io)] 

F44i  "  nr  sin2  (io]  (1  +  cos  (i0))2 

F 4 4 2  ”  T  Sin1  (io) 

F522  *  I515/32  (sin3  (iQ)  '  2  sin5  (iQ)  cos  Ci0) 

x  ? 

-  5  sin  (io)  cos“  i  *  sin  (iQ)  (-2/3 

J-  3 


+  4/3  co 5  (iQ)  +  2  cos2  (iQ) ) ) } 

F 5 4 j  =  [29.53125  sin  (iQ)  (-2-8  cos  (iQ)  +  12  cos2  (iQ) 

+  8  cos3  (i  )  -  10  cos4  (i  ))} 
o  o 

and  the  G^pq(e)  are  approximated  by  polynomials  in  e, 
depending  upon  the  range  of  value  of  e.  The  approxi¬ 
mations  used  by  NORAD  are  as  follows: 

G 2 0 - 1  =  .306  -  .44  (e  -  .64) 

(3.616  -  13 . 24  7e  +  16.29c2)  .  5<e  .  <.65 

G 2 1 1  -  (-72.099  +  331.8 19e  -  508.738e2 

+  266 . 724e3)  .65<e  <.775 

1  (-19. 302  +  117. 39e  -  228.419e2 

+  156. 591e3)  .5<e  <.65 

(-346.844  +  1582. 851e  -  2415. 925e2 

+  1246. 113e3)  . 65<e  <.775 

I  (-18.9068  +  109 . 7927e  -  214.6334e2 
G322  =  +  146 . 58 16e3)  . 5<e  <.65 

(-342.585  +  1554. 908e  -  2366. 899e2 

+  1215. 972e3)  .65<e  <.775 
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(-41.122  +  242 , 694e  -  47l,G94e 
+  313,953e3) 

(-1052.797  +  4758, 686e  -  7193. 992e2 
+  3651,9576s) 


5<e  <.65 


,65<e  <.775 


(-146.407  t  841. 88e  -  1629.0146' 


+  1083 


4356s) 


(-3581.69  +  16178. lie  -  24462. 77e' 
+  12422, 52e3) 


,  5<e  <  ,  65 


. 65<e  < . 775 


(-532.114  +  3017, 977e  -  5740.032e£ 

+  3708.276s3)  ,  5<e  <.65 

(1464.74  -  4664. 75e  +  3763. 64e2)  .65<e  <.715 

(-5149.66  +  29936. 92e  -  54087. 36e2 


31324. 566s) 


. 715<e  <.775 


(-913. 2277  +  4988. 61e  -  9064. 77e' 
+  5542. 21e3) 


. 5<e  <  .  7 


(-37995.78  +  161616. 52e  -  229838. 2e‘ 


+  109377. 94e3) 


. 7<e  <  .  78 


-V>*  XT'  %  '*%•».»«» 


17.  DISTRIBUTION  STATEMENT  (of  the  ebetrect  entered  In  Block  20,  if  different  from  Report ) 


18.  SUPPLEMENTARY  NOTES 


19  KEY  WQROS  fCondnu#  on  revets*  mtde  If  neceieery  end  tdenttly  by  block  number) 

Resonance  Tcssoral,  Sectoral  Perturbations 

Lunar  perturbations 
Solar  perturbations 
Method  of  Averaging 

20  -  ABSTRACT  f Continue  on  reverie  title  It  necemuery  end  Identity  by  block  number) 

V 

A  generalized  method  of  averaging  is  used  to  study  two  general 
problems  for  drag  free  motion  of  an  artificial  Earth  satellite. 

The  first  development  includes  perturbations  due  to  the  point 
masses  of  the  Earth,  Moon  and  Sun  as  well  as  zonal  variations. 

The  second  development  extends  the  first  to  include  tcssoral  and 
sectoral  harmonics  as  resonance  perturbations  for  two  specific 
classes  of  orbits 
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